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SUMMARY 

The critical state of vortex cores downstream of vortex breakdown has been studied. Base vortical 
flows were computed using the Reynolds-averaged, axisymmetric Navier-Stokes equations. Stan- 
dard K- e, RNG and second-order Reynolds stress models were employed. Results indicate that 
the return to supercriticality is highly dependent on the turbulence model. The K - e model pre- 
dicted a rapid return of the vortex to supercritical conditions, the location of which showed little 
sensitivity to changes in the swirl ratio. The Reynolds stress model predicted that the vortex 
remains subcritical to the end of the domain for each of the swirl ratios employed, and provided 
results in qualitative agreement with experimental work. The RNG model produced intermediate 
results, with a downstream movement in the critical location with increasing swirl. Calculations 
for which area reductions were introduced at the exit in a subcritical flow were also performed 
using the Reynolds stress model. The structure of the resulting recirculation zone was altered sig- 
nificantly. However, when area reductions were employed within supercritical flows as predicted 
using the two-equation models, no significant influence on the recirculation zone was noted. 

INTRODUCTION 

Over the past 30 years, considerable effort has been expended toward an understanding of the 
mechanisms inherent in the development and evolution of longitudinal vortices. This has been 
motivated, in part, by the desire to control and/or disable these vortices in applications such as the 
aircraft-wake-vortex hazard and submarine non-acoustic stealth. Perhaps in no application are the 
properties of swirling flows exploited to a greater extent than in the operation of gas turbine and 
industrial furnace combustion chambers. Here, a region of high swirl is induced at an entrance to 
the combustor liner, typically through a set of swirler vanes, resulting in a region of recirculating 
flow. This region acts as a fluid dynamic flameholder, providing a region of low velocity within 
which combustion may be sustained, and recirculating hot, unburned gases to the base of the 
flame. 
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This recirculating region is known as a vortex breakdown. Although several theories have 
been put forth to explain the breakdown phenomena, perhaps none are more widely recognized 
than the early works of Squire [1] and Benjamin [2]. By choosing certain functional forms for the 
base vortex flow, Squire reduced the nonlinear equations of motion (inviscid, steady) to a linear 
disturbance equation. He subsequently solved the equation to determine conditions under which a 
steady perturbation to the flow could exist. This condition, in terms of a swirl ratio, was taken as a 
limiting condition for which breakdown could occur. Benjamin examined this phenomena from a 
different perspective. He considered vortex breakdown to be a finite transition between two 
dynamically conjugate states of flow, similar to the occurrence of a hydraulic jump in open chan- 
nel flow. The two states were a subcritical state, which was defined as a flow which could support 
standing waves, and a supercritical state, unable to support standing waves. In this context, the 
work of Squire defined a critical condition marking the interface between the two states. 

The criticality condition has not received much attention from those computing numerical 
solutions to swirling flows, and flows containing vortex breakdown in particular. Of special inter- 
est is the region downstream of the breakdown. Immediately downstream, the flow is most assur- 
edly subcritical (c.f. Tsai and Widnall [3]). However, as the axial velocity recovers and the swirl 
velocity decays, the flow may return to a supercritical state at some downstream location. The 
consequences of the failure of a swirling flow to return to supercriticality has been discussed by 
Escudier and Keller [4]. In that experimental study, it was shown that the upstream influence of an 
exit contraction on vortex breakdown was substantially greater when the flow remained subcriti- 
cal compared to a flow that reverted to supercritical (upstream of the contraction). Escudier and 
Keller suggest that this phenomena might have significant consequences in the imposition of 
accurate outflow boundary conditions. 

Most swirling flows of practical interest, such as the flow within a combustor, are turbu- 
lent. For numerical calculations of turbulent swirling flows, the choice of turbulence model is of 
vital importance. It is well known that the standard K-z model does a poor job of predicting 
strongly swirling flows (c.f. Jones and Pascau [5]). One of the consequences of choosing the 
K-z model is that the wake region near the vortex centerline (downstream of breakdown) recov- 
ers much more rapidly than has been shown to occur in experiments. On the other hand, second- 
order closure models contain the physics necessary to model strongly swirling flows, and tend to 
do a better job of predicting the recovery of the axial velocity component (c.f. Jones and 
Pascau [5]). 

These varying predictive capabilities have consequences in terms of the criticality of the 
flow. That is, one would expect the K - z model to predict a return to criticality upstream of the 
position predicted by second-order Reynolds stress models. This behavior has been investigated 
to some extent by Hogg and Leschziner [6]. However, in that work swirl ratios were such that no 
recirculation zone was formed when the Reynolds stress model was employed. In addition, no 
direct calculations of the critical condition of the flow were made. However, the critical state 
(based on an inviscid analysis, c.f. Hall [7]) is not difficult to compute and thus the purpose of 
this research is to further investigate the criticality conditions of swirling flows downstream of 
turbulent vortex breakdown as predicted using several different turbulence models (K - z , RNG 
(c.f. Yakhot et al. [8]) and differential Reynolds stress). Geometries both with and without an exit 
restriction are employed. The relationship and consequences (if any) of the state of the flow (in 
terms of criticality) to the outflow restriction and turbulence model employed will be determined. 


164 



NUMERICAL PROCEDURE 


The incompressible, axisymmetric Reynolds-averaged Navier-Stokes equations are solved for the 
swirling flow within a combustor-type geometry. Although the governing equations are solved in 
general curvilinear coordinates, for purposes of brevity they are presented below in cartesian ten- 
sor form. The continuity and momentum equations are given as: 
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-r 1 « 0 
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respectively, where u- t is the mean velocity, p is the density, p is the viscosity, p is the mean pres- 
sure and t.j - u> i u ' j are ^e Reynolds stresses. 


Turbulence Models 


Although the turbulence models utilized in the present study are well documented in the litera- 
ture, the equations are included for completeness. When the K - e or RNG models are employed, 
the Boussinesq hypothesis provides an expression for the Reynolds stresses is terms of the gradi- 
ents of the mean flow as: 
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where v, is the turbulent viscosity, K is the turbulent kinetic energy and S- is the strain rate. The 
turbulent viscosity is expressed in terms of K and the dissipation rate e as: 



(4) 


Transport equations for K and e, respectively, are written as: 
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In the case of the standard K-e equation R - 0. For the RNG model (see Yakhot et al. [8]): 
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where r) - SK / e and S = ( 25^5 ») 1/2 . 

It remains to specify the constants in the above equations. For the K - e model the stan- 
dard values for boundary layer flows (C^ - 0.09, C ei - 1.44, C t2 - 1.92, - 1.0 and 
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a - 1.3) have been taken. For the RNG model, theoretical analysis yields that C e2 - 1.68, 
C tl - 1.42, a k - o £ - 0.72, ri 0 - 4.38 and p - 0.012. 

The Reynolds stress model involves the solution of transport equations for the individual 
Reynolds stresses The following equations, employing the closure assumptions of Gibson 
and Launder [9], and Launder [10], are solved: 
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The production term is computed as: 
du ; _ 


p v - 
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The pressure/strain and dissipation terms are modeled as: 
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with values C 3 - 1.8 and C 4 - 0.06 assigned to the constants. 

For each model wall functions based upon the assumption of a fully developed equilib- 
rium turbulent boundary layer are utilized in the near wall region. This approach is deemed suit 
able since the physics of the problem are not dominated by near wall phenomena. 


Solution Procedure 

The above equations were solved using the commercial code FLUENT [11]. FLUENT uti- 
lizes a pressure-based control volume technique. Second-order upwind interpolation is used to 
provide values of variables on cell faces. Pressure-velocity coupling is implemented using the 
SIMPLEC algorithm [12]. Convergence of the solution is assumed when the sum of the normal- 
ized residuals for the conservation equations is decreased to a minimum of 1.0 x 10 -3 . (The resid- 
ual for a given equation consists of the summation of the unbalance in the equation for each cell in 
the domain.) Comparisons revealed that solutions converged by an additional factor of two were 
virtually indistinguishable. Since the above techniques are well known and widely discussed in 
the literature, they will not be elaborated upon here. 

Geometry and Boundary Conditions 

The geometric configuration is that of a prototypical combustor. The geometry near the expan- 
sion/breakdown region is shown with the computational grid superimposed in Figure 1. The 
domain was discretized using 170 (axial) and 40 (radial) grid points. Grid points were clustered 
near the breakdown region. The length of the domain is 40/i, where h is the radius at inflow. The 
expansion of area ratio 4:1 takes place over a length of 3 h. Note that the i=constant lines in the 
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physical grid occur at constant x — this is to facilitate calculation of the flow criticality. Addi- 
tional calculations were made utilizing a 114 x 27 grid. No significant changes in the solutions 
were noted. 

At inflow, a solid body rotation of the form w - Qr has been specified (where w repre- 
sents the circumferential velocity). The value of Q was assigned values ranging form 0.75 to 1.5. 
The axial velocity has been taken as uniform, that is u - 1.0, and the turbulence intensity was set 
at 10%. At the outflow boundary, zero streamwise gradient conditions were enforced. The appro- 
priateness of these conditions is confirmed by examining the distribution of the velocity contours 
(as shown in the Results section) near the outflow boundary. In addition, for each calculation the 
Reynolds number, based on the axial velocity and duct radius at inflow, was 100,000. 


Criticality Calculations 


Determination of the criticality of the flow is based the solution of the 


ential equation (c.f. Benjamin [2] or Hall [7]): 
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where k - rw is the circulation and F c is a quasi-cylindrical perturbation shape function. 
Assumptions inherent in the above equation are that the fluid is inviscid, and that the flow is 
steady and axisymmetric. Thus, the criticality condition, as defined by Benjamin or Hall, is con- 
cerned only with the propagation or existence of axisymmetric waves on inviscid cores. However, 
we shall utilize the theory to predict the ability of high Reynolds number turbulent cores to sup- 
port axisymmetric waves. __ 

The radial distributions for u and K are available at any axial location from the mean flow 
<g^culations. The above equation is solved subject to the boundary conditions F c - 0 and 

c « constant at r - 0. The flow is subcritical if the solution curve passes through zero in the 

interval 0<r<2h, supercritical if the solution does not pass through zero, and critical if the solu- 
tion is zero at both r - 0 and r - 2h (where 2h is the duct radius). The above equation is 
solved utilizing a second-order accurate modified Euler technique. 

RESULTS 


The behavior of each of the turbulence models in predicting the location of the return to criticality 
of the vortex is presented first. Following this, contours of constant axial velocity are examined, 
and the effects of outlet restrictions on the flow are discussed. 

Mean flow calculations were performed for Q=0.75, 1.125 and 1.5, utilizing the K- e, 
RNG and Reynolds stress turbulence models. The axial location at which the vortex returned to a 
supercritical state was computed by solving an ordinary differential equation, as described earlier. 
In each case, the vortex was supercritical upstream and subcritical immediately downstream of 
the breakdown location. The results for return to supercriticality are summarized in Figure 2. The 
K-z model predicted a return to supercritical conditions at locations considerably upstream of 
that predicted by the other models. In addition, and contrary to what one would expect, th tK-z 
model predicted that the location at which the vortex returned to a supercritical state was not 
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affected by changes in the swirl level at inflow. Conversely, the RNG model showed an expected 
sensitivity to increases in swirl. That is, as the swirl level at inflow was increased, the critical 
location moved downstream. In the case of the Reynolds stress model, the vortex remained sub- 
critical to the exit for ail swirl levels. We note that Q - 0.75 represented the approximate mini- 
mum swirl level for which a recirculation zone was formed using the Reynolds stress model. 
Consequently, calculations of vortex breakdown flows using Reynolds stress models will likely 
involve outflow boundary conditions imposed on subcritical flows (due to practical limitations on 
the length of the computational domain). However this need not be the case for the two-equation 
models. 

Escudier and Keller [4] have shown that the shape and internal structure of the recircula- 
tion region is strongly influenced by outlet restrictions if the restriction is imposed within a sub- 
critical flow. To examine whether or not these effects are observed in numerical calculations, 
results have been computed using the Reynolds stress model for area reductions of 19% and 36% 
with Q - 1.5. The reductions extended over 33.67 h 40 h. Results are presented in terms of 

axial velocity contours in Figures 3a-c for the unrestricted case, the 19% reduction, and the 36% 
reduction, respectively. (For purposes of clarity, the contour plots have been scaled by a factor of 
3 in the radial direction. In addition, the inner-most contour levels within the bubble structure rep- 
resent the level u - -0.1.) It is clear that the 36% area reduction has a large effect on the shape 
of the aft portion of the bubble near the vortex centerline. The zero axial velocity contour under- 
cuts this portion of the bubble, in effect lifting the recirculation zone off the axis. Note however, 
that the forward portion of the bubble remains virtually unaffected. In addition, a strong jet-like 
vortex core exists downstream of the breakdown. These features are in very good agreement with 
those observed by Escudier and Keller [4] in their experimental work. For the 19% area reduction 
case, the primary effects concern the rate of recovery of the axial velocity downstream of the 
breakdown. The shape and internal structure of the bubble remains similar to that of the unre- 
stricted case. This result is also consistent with the results of Escudier and Keller. We do note that 
the geometry of Escudier and Keller was somewhat different than that employed in this study. In 
their study, a solid inner cylinder was included at the inflow plane, and a step expansion rather 
than a gradual expansion was used. Numerical convergence problems prevented our use of that 
geometry. However, the structure of the recirculation region for confined flows undergoing vortex 
breakdown does not appear to be overly sensitive to the inlet geometry. 

It is of interest that the downstream measurements by Escudier and Keller were 0.39L 
upstream of the area reduction (where L was the distance between the inlet and the reduction). 
Thus, it is quite possible that for their cases in which the area reductions were very large, the flow 
actually returned to a supercritical state upstream of the restriction (but downstream of their last 
data point). That is, due to continuity the mean axial flow velocity increases with the square of the 
area reduction; however conservation of angular momentum dictates that the swirl velocity 
increase in a linear manner. Thus, the ratio of swirl to axial velocity generally decreases, resulting 
in a possible return to supercritical conditions. In fact, the authors found that for the low swirl 
case (£2 - 0.75) area reductions on the order of 20% did result in the flow returning to a super- 
critical state slightly upstream of the area reduction. However, for the high swirl cases shown in 
Figures 3b-c, the vortex remained in a subcritical state to the exit. 

Contours of constant axial velocity for the K - t and RNG models (without exit restric- 
tions) are shown in Figures 4a-b, respectively for the case Q - 1.5. The differences in the axial 
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velocity distributions downstream of the breakdown region as predicted by these models and the 
Reynolds stress model are quite large. The rapid increases in the axial velocity for the K - e and 
RNG models (accompanied by equally rapid decreases in the swirl velocity) account for the rapid 
return to supercriticality. In addition, the internal structures of the bubbles predicted by these two- 
equation models differ considerably amongst themselves, and with that produced by the Reynolds 
stress model. Calculations were also made for the K - e and RNG models 36% area reductions. 
For the sake of brevity, we show only the results from the RNG model in Figure 4c. (Recall that 
for the two-equation models the vortex was supercritical upstream of the restriction.) As the fig- 
ure reveals, virtually no difference in the axial extent or shape of the recirculation region resulted 
(when compared with the case in which no restriction was employed). This further confirms the 
experimental results of Escudier and Keller, and highlights the extreme sensitivity of strongly 
swirling flows to the turbulence model employed. 

CONCLUSIONS 

It is clear from the results of this study that wide differences exist in the predictions of two-equa- 
tion and Reynolds stress turbulence models for return to supercriticality of swirling flows down- 
stream of vortex breakdown. The results of the Reynolds stress model are in better agreement 
with experimental results for swirling flows (in similar geometries) as described by Escudier and 
Keller [4] than results predicted using two-equation models. Thus, the study further highlights the 
inappropriateness of using two-equation models to predict strongly swirling flows. The sugges- 
tion by Escudier and Keller that outflow restrictions might have a drastic effect on the structure of 
the breakdown as predicted through numerical solutions to the equations of motion was con- 
firmed. Results revealed that for relatively large area reductions (on the order of 36%), the struc- 
ture of the recirculation region may be greatly affected. However, for lesser reductions, the shape 
of the recirculation zone may be very similar to that resulting from the unrestricted geometry. 
Consequently, it does not appear that the requirements for the specification of outflow boundary 
conditions need be significantly more stringent for subcritical swirling flows than for supercritical 
flows, or for flows without swirl. This is fortunate — the persistent nature of the subcritical flow as 
revealed by the experiments of Escudier and Keller [4], and as predicted by the Reynolds stress 
model would severely restrict the predictive capability of many engineering-type calculations. 

The authors plan future work in the area of combusting flows. For these calculations, 
accelerations in the axial velocity due to decreases in density should considerably alter the critical 
state of the flow. 
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3. Contours of constant axial velocity; Reynolds stress model, Q - 1.5; contour levels from • 0.1 
to 1.0 in intervals of 0.1 (geometry scaled by a factor of 3 in the radial direction). 

a) No outlet restriction 

b) 19% area reduction 


c) 36% area reduction 
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4. Contours of constant axial velocity; Q - 1.5 (geometry scaled by a factor of 3 in the radial 
direction). 

a ) K-t model; contour levels from -0.2 to 1.0 in intervals of 0.1. 

P b) RNG model; contour levels from -0.4 to 1.0 in intervals of 0.1. 

c) RNG model; 36% area reduction; contour levels from -0.4 to 1.0 in intervals of 0.1. 


i 


173 


